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ABSTRACT 

Long-exposure Suzaku X-ray observations of the nearby active galaxy NGC 405 1 from 2005 
and 2008 are analysed, in an attempt to reach a self-consistent understanding of both the spec- 
tral variability on long timescales and the broad-band variability at high time resolution. The 
techniques of principal components analysis and a maximum likelihood method of power 
spectrum analysis are used. In common with other type I active galactic nuclei (AGN), the 
spectral variability is dominated by a varying-normalisation power-law component together 
with a quasi-steady, hard-spectrum offset component that contains FeK atomic features. 
NGC 4051 displays a strong excess over a power-law at energies above 20keV, some frac- 
tion of which also appears to vary with the power-law continuum. The fast timescale power 
spectrum has a shape generally consistent with previous determinations, with the previously- 
known dependence on broad-band photon energy, but in the new data significant differences 
are found between the low and high flux states of the source, demonstrating the power spec- 
trum is non-stationary. Frequency-dependent time lags between the hard and soft bands of up 
to 970 ± 225 s are measured. The existence of the observed time lags excludes the possibility 
that the hard spectral component originates as reflection from the inner accretion disk. We in- 
stead show that the time lags and their frequency- and energy-dependence may be explained 
simply by the effects of reverberation in the hard band, caused by reflection from a thick 
shell of material with maximum lags of about 10 000 s. If the reflecting material surrounds the 
AGN, it extends to a distance about 1.5 x 10 14 cm, 600 gravitational radii, from the illuminat- 
ing source and the global covering factor is C g > 0.44, confirming previous suggestions that 
type I AGN have high covering factors of absorbing and reflecting material. Given the spectral 
and timing similarities with other type I AGN, we infer that this source structure is common 
in the type I population. 

Key words: accretion, accretion discs - galaxies: active - X-rays: galaxies - X-rays: individ- 
ual: NGC 4051 



1 INTRODUCTION 

Given the impossibility of imaging the innermost regions of ac- 
tive galactic nuclei (AGN) at X-ray energies, the principal tools 
for studying accretion in this waveband have been either modelling 
of X-ray spectra or studying the rapid X-ray variability of nearby, 
bright AGN. Owing to signal-to-noise limitations spectral analy- 
sis has tended to c oncentrate on fitting simple composite models to 
mean spectra (e.g. lMiniutti et alJl2007D . For the same reason, tim- 
ing analyses have studied the variability on short timescales but in 



broad bandpasses of photon energy (e.g. lMcHardvrtlii]|2004l) . To 
date, conclusions from these studies have not been combined in 
a self-consistent way. However, there is substantially more infor- 
mation available in the spectral variability than is contained in the 
mean spectrum alone, and now that some long-exposure observa- 
tions have been made, those data may be time-sliced and the spec- 
tral variability used as a powerful constraint on the models (e.g. 
iMiller et alj|200l IMiller. Turner & Reevesll2008l) . Long exposure 
observations have also revealed interesting variability signatures on 
short timescales (see below). In this paper we aim to carry out both 
spectral variability and timing analysis on long-duration Suzaku ob- 
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Table 1. Summary of the three Suzaku observations, giving the ID number, dates, the detectors used in the analysis, total duration, XIS on-source time, the 
mean XIS count rate in the 0.3-10keV range, and the mean flux in the 2-10keV range. The 2005 count rate has been scaled to the value for two detectors, for 
consistency with the 2008 observations. 

ID dates detectors duration on-time mean count rate mean 2-10 keV flux 

/s /s /s _1 /10" 11 ergs~ 1 cm -2 

700004010 2005 Nov 10-13 XIS 0,2,3, PIN 222592 119618 (L93 (187 

703023010 2008 Nov 6-12 XIS0.3, PIN 494834 274531 4.19 2.41 

703023020 2008 Nov 23-25 XIS0.3, PIN 161695 78394 2.87 1.79 



servations of the nearby AGN NGC 4051, and to reach a consistent 
understanding of the source behaviour. 

There are two phenomena of particular interest in each of the 
"spectral variability" and the "timing" domains. The first is that it 
is now well-established that type I AGN frequently show a broad 
"red wing" of emission extending below the 6.4 keV Fe Kq emis- 
sion line that appears as an ex cess above a model power-law contin- 
uum (e.g. Tana ka et al ■I 19951) . and this as been widely interpreted as 
being Fe Ka emission from within a few Schwarzschild radii of the 
black hole. Model fits to mean spectra infer a redshift so large that 
the emission would need to come from within the innermost stable 
circular orbit of a Schwarzschild black hole, lea ding to claims of th e 
detection of black hole spin in these AGN (e.g. lWilms et al.ll200ll) . 
However, a particular problem is that any Fe Kq emission should 
vary in phase with the illuminating continuum, whereas the red 
wing appears to be largely invariant while th e illuminating source 
appears to exhibit large time variations (e.g. Ilwasawa et al]| 19961 : 
Vaughan & Fabian 2004). The problem is now known to extend to 
high energies: above 20 keV many AGN have an excess of contin- 
uum compared with extrapolation of a 2-10keV power-law model, 
and this has been supposed to be continuum reflection (Compton 
scattering modified by the photoelectric opacity of the gas) from 
the same material that produces the red wing. Indeed, this "hard ex- 
cess" is likewise known to exhibit little variability compared with 
the 2-10keV band, p ointing to a common origin with the red wing 
dMiniutti et alj|2007T) . However, its flux is so large that it cannot 
readily be ascribed to reflection unless the illum inating source is 
signifi cantly obscured along our line of sight (e.g. iTerashima et aU 
120091) . 

To overcome thes e problems. iFabian & Vaughanl d2003h . 
iMiniutti etaTI d2003l) and lMiniutti & Fabianl d2004l) haveproposed 
a "light-bending" model in which the effect of light following 
geodesies near the black hole mean that a reflection component 
from the inner accretion disk may be made relatively invariant, if 
it is illuminated by a source whose position changes with respect 
to the disk. The model requires both that the illuminating source 
is compact, comparable in size to the event horizon, and moves 
towards and away from the disk. As it moves, the flux received 
by a distant observer varies as light is removed from the line of 
sight. In this model, it is the change in position of the primary 
source that dominates the observed variability, rather than intrin- 
sic fluctuations. The light bending could also explain the h i gh ap- 
par ent reflectivity in the hard band . However. iMiller et alj ( 2008) 
and lMiller. Turner & Reeve's] d2009l) have shown that the red wing 
and the hard excess in MCG-6-30-15 may instead be explained by 
the spectrum expected from clumpy, "partial-covering" absorbers, 
and that a model in which the absorber covering fraction varies pro- 
vides a good description of the observed spectral variability, as well 
as of the mean spectrum of this AGN. Both Compton-scattering re- 
flection and absorption shape this hard spectral component, but the 
relative importance of scat tering/reflection a nd absorbed transmis- 
sion is not yet established (Miller et al. 200^)- 



The second important phenomenon occurs in the timing do- 
main. There are now a significant number of AGN where time lags 
in the range 10s to 1000s of seconds between different spectral 
bands are seen, in the sense that hard - band photons lag soft-band 
1 Papadakis. Nandra & Kazanasl 200 ll Vaughan. Fabian & Nandral 



20031IM cHardy et alj|2004l. lMarkowitz!l2005L lArevalo et alj|2006l. 
Markowitz et al.ll2007l lArevalo e t al. 2008). Importantly, the lags 
are dependent on the frequency of the source variation, in the sense 
that if a source's variations are decomposed into Fourier modes, 
the lag between hard and soft energy bands increases with the pe- 
riod of those modes, and also increases with the separation in en- 
ergy of the bands. The lags and their frequency-dependence have 
been explained by a model comprising perturbations propagating 
inwards on an accretio n disk where harder X-ra y emission is emit- 
ted from smaller radii ( Arevalo & Uttley 2006). Although the ini- 
tial observations were consistent with a linear relationship with 
Fourier period, recently, more complex behaviour has been seen, 
leadin g to a hypothesis that fluctuations arise in multiple compo- 
nents (McHardy et al. 2007). However, this explanation of rapid 
fluctuations propagating through the accretion disk is not consis- 
tent with the light-bending model, which has the primary varia- 
tions caused by motion of the hypothesised compact source above 
the accretion disk. It is our aim in this paper to reconcile the phe- 
nomena of spectral variability and time lags in a detailed study of 
NGC 4051. 

NGC 4051 is a well-studied, nearby narrow-line Type I AGN 
with redshift z = 0.0023. If the redshift were purely due to cos- 
mological Hubble flow its distance would be 9.3 Mpc for Hubble's 
constant Ho = 74kms _1 Mpc -1 , but for such a nearby galaxy 
the Tully- Fisher distanc e is probably more reliable, placing it at 
15.2 Mpc dRusselll2004h . It has a black hole mass determined from 
optical reverberation mapping of Mbh = l-7lo's2 x 10 6 M0 
dDennev et alj 120091) . It has long been known to be highly vari- 
able on short timescales at X-ray energies, with spectral vari- 
ability that is correlated with flux { Lawrence et al. 1985[ 1 19871 
i Matsuoka et"ai1ll990l. fkunieda et alj|l992l iMcHardv et aljfl995l 
lLamer et alj 120031) . The timing properties of X-ray observations 
with the Rossi X-ray Timing E xplorer (RXTE) and X MM-Newton 
have been studied in detail by [McHardy et al. ( 2004) . The X-ray 
spect rum has previously been studied by iGuainazzi et aL I d 19961 
1998)and |pounds et al J d2004l) . among others. Of particular interest 
is a 2005 Suzaku observation obtained when the source was in a low 
state that reveals a str ong excess of emission above a power-law at 
energies E > 20keV dTerashima etai]|2009l) . Overall, the consen- 
sus from recent analyses is that NGC 405 1 has low states where the 
primary power-law source has largely disappeared and where the 
remaining emission appears to be dominated by reflected emissio n 
dGuainazzi et aljll998|: Pounds et alj|2004l : ITerashima et~ai]|2009h . 
ITerashima et alj (2009) found that partially-covered reflection was 
required to fit the hard excess, rather than supposing this to be 
inner-disk gravitationally-redshifted reflection. The high variabil- 
ity on short timescales means that this AGN provides a powerful 
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testbed for using spectral variability and power spectrum analysis 
methods to untangle the source emission components. 

We now consider both the 2005 Suzaku observation and two 
new Suzaku observations obtained in 2008, and analyse both the 
time-dependent spectrum and the high time resolution behaviour. 



2 OBSERVATIONS AND DATA REDUCTION 

The Suzaku data analysed here were reduced as described by 
iTurner et alJfeOlti) and Lobban et al. (in preparation). The Suzaku 
X-ray Imaging Spectro meter (XIS, Kovama et al. 2007) comprises 
four X-ray telescopes dMitsuda et aTl r2007) each with a CCD de- 
tector. XIS CCDs 0, 2, 3 are configured to be front-illuminated and 
yield useful data over ~ 0.6 — 10.0 keV with energy resolution 
FWHM~ 150 eV at 6keV. XIS 1 is a back-illuminated CCD and 
has an enhanced soft-band response down to 0.2 keV but lower area 
at 6 keV than the front-illuminated CCDs, as well as a higher back- 
ground level at high energies. XIS 1 was not used in our analysis. 
Use of XIS 2 was discontinued after a charge leak was discovered 
in Nov 2006, so we used data from XIS 0, 2 and 3 in 2005, and data 
from XIS and 3 in 2008. 

Our analysis used three observations from 2005 Nov 10-13, 
2008 Nov 6-12 and 2008 Nov 23-25. The data were reduced using 
v6.4. 1 of HE ASOFT and screened to exclude periods, i) during and 
within 500 seconds of the South Atlantic Anomaly, ii) with an Earth 
elevation angle less than 10° and iii) with cut-off rigidity > 6 GeV. 
The source was observed at the nominal pointing for the XIS. The 
front-illuminated CCDs were in 3 x 3 and 5x5 edit-modes, with 
normal clocking mode. For the XIS we selected events with grades 
0,2,3,4, and 6 and removed hot and flickering pixels using the SIS- 
CLEAN script. The spaced-row charge injection was used. The XIS 
products were extracted from circular regions of 2.9'radius with 
background spectra from a region of the same size, offset from the 
source (avoiding the calibration sources at the edges of the chips). 
The response and ancillary response files were created using XIS- 
RMFGEN V2007 MAY and XISSIMARFGEN V2008 MAR. 

Suzak u also has a non-imagi ng, collimated Hard X-ray Detec- 
tor (HXD, lTakahashietal]|2007l) whose PIN instrument provides 
useful data over 15-70 keV for bright AGN. NGC 4051 is too faint 
to be detected by the HXD GSO instrument, but was detected 
in the PIN. For the an alysis we used the model "D" background 
( Fuka zawaetaiT2 009 ) . As the PIN background rate is strongly vari- 
able around the orbit, we first selected source data to discard events 
within 500s of a South Atlantic Anomaly passage, we also rejected 
events with day/night elevation angles > 5°. The time filter result- 
ing from the screening was then applied to the background events 
model file to give PIN model background data for the same time in- 
tervals covered by the on-source data. As the background events file 
was generated using ten times the actual background count rate, an 
adjustment to the background spectrum was applied to account for 
this factor. HXDDTCOR V2007 MAY was run to apply the deadtime 
correction to the source spectrum. As the PIN deadtime correction 
varies, for the time-dependent analysis below, the correction was 
calculated at 128 s time resolution and the mean correction in each 
analysed time interval was applied to the source counts. 

The observations are summarised in TableQ] which also gives 
the XIS count rates and the flux integrated over 2-10keV. The 2005 
observation caught NGC 4051 in a historical low state, while in 
2008 the source flux was more typical, and close to the valu e in the 
2001 XMM-Newton data analysed bv lMcHardv et alj J2004I) . 
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Figure 1. Combined XIS and PIN spectra of NGC 4051 from 2005 (lower 
curve) and 2008 Nov 6-12 (upper curve), plotted as Ef(E), unfolded 
against a power-law with photon index T = 2. XIS data have been binned 
at the instrumental HWHM, PIN data at intervals A log 10 E = 0.1. Photon 
shot noise errors are shown. No cross-calibration correction between XIS 
and PIN instruments has been applied in this figure. 



3 SPECTRAL VARIABILITY ON TIMESCALES 
LONGER THAN 20 KS 

3.1 Principal components analysis 

Our first step to understanding the variations in NGC 4051 is to 
carry out a principal components analysis (PCA) of the spectral 
variations. This is a method of decomposing the variations into or- 
thogonal components, eigenvectors. By ordering the eigenvectors 
according to their eigenvalues we can extract the principal modes of 
variation. It is particularly well suited to variations that are additive 
in nature (such as variations in partial covering fraction by a clumpy 
absorb er). We use the method described in detail by Mill er et al.l 
(2007), in which spectral bins are created of width equal to half the 
energy-dependent instrument resolution. Spectra are averaged over 
and sampled at time intervals that are sufficiently long that varia- 
tions between time samplings are dominated by intrinsic variations 
rather than shot noise. The method uses singular value decompo- 
sition to cope with the problem that there are more spectral points 
than there are available time bins: in this case the leading eigenvec- 
tors are still defined and may be extracted. 

In interpreting the results, there should be no expectation that 
the orthogonal eigenspectra produced correspond uniquely to any 
physical component in the AGN. However, what has been f ound in 
previous analyses of other AGN dMiller et al1l2007l, 120081) is that 
only a small number of variable components are required to de- 
scribe the spectral variability, which is indicative that additive pro- 
cesses are at work. In those previous analyses we concentrated on 
eigenvector one, which describes the primary source variation, and 
subsumed the rest of the source spectrum i nto an "offset" com po- 
nent which is esse ntially steady in time (see lMiller et alj|2007l and 
iTurner et alj2007l for more detailed discussion). The offset compo- 
nent is not uniquely defined, because we are free to add or subtract 
arbitrary amounts of eigenvector one without changing the decom- 
position. However, we may define a minimum offset component by 
requiring that no spectral bin have negative flux and a maximum 
offset component by requiring that no spectral bin be higher than 
the lowest observed value. 
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Figure 2. Principal components analysis of the combined 2005-8 dataset for NGC4051: (left) the offset component, showing the possible range of this 
component (upper and lower curves) and the component of cosmic X-ray background that has been subtracted (solid curve). For clarity, errors are shown only 
every fifth point on one of the offset spectra; (right) eigenvector one. 



3.2 Spectral variability 

Fig-U shows the mean spectra from 2005 and 2008 Nov 6- 12, plot- 
ted as unfolded spectra against a power-law with photon index 
r = 2. In the figure, a cosmic background model has been sub- 
tracted from the PIN data by using XSPEC v 11.3.2ag to generate 
a spe ctrum from a cosmic X-ray background model dGruber et al.l 
1999) normalized to the 34' x 34' PIN field of view, which was then 
combined with the PIN instrument background file. No correction 
for the calibration difference between the XIS and PIN instruments 
has been applied in Fig.[T] this has the effect of making the PIN 
data too high by 16 percent relative to the XIS data. The basic spec- 
tral variability mode of the source is clearly seen, with the hard, 
low-state and higher softer spectrum at the two epochs. There is 
a significant hard excess above 20keV, which has varied between 
the epochs. In 2008 there is also a "red wing" as an excess above 
a power-law continuum below the 6.4keVFeKa line: in the 2005 
data this can be seen to have developed into a rather different con- 
tinuum shape. Previous analyses have assumed the red wing to be 
relativistically broadened Fe Ka emission (e.g. lLamer et alj|2003l) 
but this interpretation is not obviously supported by the PCA. 

The PCA of the observations in the energy range 0.4-50 keV is 
shown in Fig.|2] The time sampling was 20 ks: although this source 
is highly variable on shorter timescales, it is important to avoid the 
regime where counts in each spectral bin are dominated by shot 
noise. The full 2005-8 dataset was used, incorporating both XIS 
and PIN data. XIS data were binned at the spectral resolution half- 
width at half-maximum. To achieve adequate signal-to-noise in the 
PIN data, coarser spectral binning was used (Alog 10 E = 0.1). 
The minimum and maximum limits on the offset component are 
shown, together with the contribution of the cosmic X-ray back- 
ground that has been subtracted from the PIN offset component. 
The PIN data have been reduced by the 16 percent cross-calibration 
factor between the XIS and PIN instruments. Uncertainties on each 
sp ectral point were ob tained by a Monte Carlo method as described 
by Mil ler et all |2007). Eigenvector one accounts for 81 percent of 
the variance in this energy band. Eigenvector two accounts for a 
further 1 1 percent of the variance, and although noisy this compo- 
nent appears dominated by a residual anticorrelation between the 
PIN and XIS data, which when considered with the eigenvector one 
spectrum implies there is a minor element of uncorrelated variation 
between the XIS and PIN bands. In this paper we shall concentrate 



on the gross spectral variability characterised by eigenvector one 
and the offset component. 

The PCA shows the same generic features as the other type I 
AGN that it has been applied to (e .g. Mrk 766. IMiller et al.|[2007l 
MCG-6-30-15, IMiller et al.l 120081) . The hard offset component 
dominates the PIN band, and contains atomic signatures, partic- 
ularly of the FeK edge and FeKa emission line, indicative of 
the importance of atomic processes in shaping this component. 
Eigenvector-one is essentially a power-law, but with some "warm 
absorber" signatures imprinted on it in the soft band, and, unusu- 
ally, evidence of both an excess variable component above 20 keV 
and weak evidence for broadened Fe Ka emission, with an equiva- 
lent width of about 66I37 eV. The width of the broad component is 
not well constrained, but there is no evidence for any excess vari- 
able component below 6 keV, and in common with other AGN, the 
"red wing" component is described by the constant offset spectrum: 
the red wing feature in the 2008 data is implied by the PCA to be 
simply a sum of the very hard-spectrum offset component and the 
power-law-like eigenvector one. 

The interpretation of eigenvector-one is fairly straightforward, 
as it implies that the primary continuum has a power-law form 
of invariant index, whose amplitude varies. The amplitude varia- 
tions could either be intrinsic to the source producing the power- 
law, or could arise from changes in covering fraction of an exter- 
nal absorber. The variable hard excess could either be flux trans- 
mitted through a partially-c overing, high colum n density absorber, 
such as found in PDS456 l iReeves et al]|2009h and 1H 0419-577 
dTumer et al.ll2009h . or it could arise as reflection from optically- 
thick, or nearly so, material. If there are high column density 
partial-covering absorbers, we expect to see b oth transmitted flu x 
at high energy and Compton-scattered flux (see IMiller et al.l2009h . 

The offset component, which essentially describes the source 
in its lowest observed flux state, may be explained either as be- 
ing dominated by reflection from optically-thick material, or as 
being the residual signal that arises from par tial-covering varia- 
tion s. These poss i bilitie s are discussed more bv lMiller et afl l l2007l) 
and lTurner et alj d2007h . where the point is made that it is funda- 
mentally not possible to distinguish between these possibilities by 
continuum fitting alone. Although spectral models have previously 
tended to be characterised as being dominated by either one or other 
of absorption or reflection, it now seems likely that both physical 
processes make significant contributions: absorbing zones of high 
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Figure 3. XIS light curves, 0.3-10keV, for the three Suzaku observations, sampled at 256 s: (left) 2005, (centre) 2008 Nov 6-12, (right) 2008 Nov 23-25. The 
XIS 0,2,3 observations of 2005 have been scaled to the flux expected in two detectors, to allow comparison with the XIS 0,3 measurements in 2008. Gaps 
indicate periods of no useful data. Time values are set to zero at the first period of useful data in each observation. 



column de nsity also contribu te a significant flux of Compton scat- 
tered light (Miller etalj2009h . We return to discussion of the inter- 
pretation of the spectral variability in section[5] after consideration 
of the variability on shorter timescales. 

Finally, the FeKa emission is particularly interesting The 
narrow 6.4 keV line, produced by Fel-XVII, appears only on the 
offset component, and hence must have had an invariant flux in the 
2005 and 2008 epochs. The implication is that the gas responsible 
for the line experienced the same ionising continuum in 2005 as in 
2008, despite the observed 7-10 keV flux being a factor 2.5 higher 
in 2008 than in 2005. We discuss later the interpretation that the dip 
into the 2005 low state was caused by obscuration along our line of 
sight. 



3.3 Summary of key spectral variability results 

We summarise the key points from the above, to aid the later dis- 
cussion. 

(i) NGC 405 1 shows the same general form of spectral variabil- 
ity as other type I AGN, with a soft varying power-law and a less 
variable hard spectrum component that creates a red wing below 

6.4 keV and PIN-band excess above the power-law in the mean 
spectrum. The 2005 low state is dominated by this h ard compo- 
nent. T he new analysis confirms the conclusions of Poun ds et al.l 
d2004l) and lTerashima et al.l fe009l) . 

(ii) The 6.4 keV Fe Ka narrow emission line was unchanged be- 
tween 2005 and 2008 d espite a significant c hange in continuum 
flux, as also discussed bv lPounds et alj d2004l) . 

(iii) There is also a component of variable hard excess at 
E > 20keV that varies with the power-law continuum on 20 ks 
timescales, with evidence for an associated Fe Ka line component 
that is moderately broadened (i.e. not relativistically broadened), 
by either Doppler or Compton scattering effects. 



1 iTurner et al.1 d2010h also discuss the companion line visible at 5.44 keV 
in the offset spectrum. 



4 POWER-SPECTRUM AND TIME DELAY 
MEASUREMENT 

4.1 Maximum-likelihood estimation of the power spectrum 
and time delays 

We now wish to complement the spectral variability analysis with 
analysis of the spectral variations in broader energy bands on 
shorter timescales. In the case of Suzaku obs ervations it is not possi- 
ble to use periodogram-type methods (e.g. Papadakis & Lawrence! 
1 19931) because of the substantial gaps in the data caused by the 
filtering of events and the short orbital timescale. TableQ] shows 
that the fraction of usable data is in the range 0.48-0.55. Further- 
more, the orbital period imposes regular gaps of period 5753 s on 
the dat a. Instead we use a maximu m likelihood method, based on 
that o f Bond , Jaffe & Knoxl dl998l) . which is described more fully 
in the context of X-ray time series analysis by Miller (in prepara- 
tion). The reader is referred to that paper for a full description of the 
method, but briefly, we start by defining the power spectral density 
(PSD), P{v), as the Fourier transform of the autocorrelation func- 
tion, A(t), 



A(r) 



&vP(v) cos(27t^t) sine' 



\2u N J 



(1) 



where is the sampling Nyquist frequency (we assume uni- 
form sampling in time) and the sine 2 term arises from the effect 
of binning the c ounts in bins of duration equal to the time sam- 
pling (see also lEfstathiou & Moodvl 1200 ll) . We do not construct 
a periodogram, but instead use an iterative maximum likelihood 
method to find the best-fitting coefficients for a model PSD given 
the measured time series. We adopt a gaussian likelihood function, 
£, which thus takes the form 



C = 



(2it) N / 2 \C\ 1/2 



exp 



(2) 



for a time series containing N data samples, where C is a N x 
N model covariance matrix and A is the time series of fractional 
fluctuations of length N. To find the best-fitting PSD we compute 
the expected covariance matrix for an initial set of model parameter 
values using equationQ] and iterate. Because we can calculate the 
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derivatives of the likelihood function with respect to the model, 
fro m equationQ] we can use a Newton-type method as described 
bv lBond et al.l d!998h . and the iterations typically converge within 
5 iterations to a stable solution irrespective of the number of model 

parameters. 

The PSD model we adopt follows iBond et alj Jl998h in defin- 
ing a set of "bandpowers", which we choose to be stepwise in vP v : 

A(t) = 2 ^ (uiPi) cos(2ttvt) sine 2 (3) 

where each bandpower has PSD amplitude Pi in the frequency 
range fi,i-i>2,j- Within each band the PSD is constant in vP(y), 
which helps reduce sharp steps in the PSD at bandpower edges 
compared with defining bandpowers uniform in P(y). As part of 
the likelihood maximisation procedure, we calculate the Fisher ma- 
trix for the PSD coefficients, allowing us to have estimates of the 
errors and their covariance. Those errors include the sampling un- 
certainty that arises because we only observe a finite time series. 
The likelihood estimation is in the time domain, so the errors on 
the PSD coefficients are well-understood even in the presence of 
a complex sampling window function. The method is immune to 
gaps in the data, periodic or otherwise (if the time sampling is such 
that one particular Fourier mode is badly sampled, this is simply 
reflected in the corresponding uncertainty). Multiple datasets span- 
ning arbitrarily long time baselines may be trivially combined, so 
we are able to analyse the three Suzaku observations as either indi- 
vidual or combined time series. 

The likelihood function that is maximised includes terms in 
the autocorrelation function for both modes intrinsic to the source, 
described by P(v), and Poisson shot noise. Thus the likelihood 
method automatically corrects for Poisson noise in the data and 
is accurate even in the low count rate regime. Because the shot 
noise is put into the likelihood model, the uncertainties on the PSD 
that we measure include the statistical uncertainty arising from shot 
noise. Although correction for high-frequency aliasing may also be 
included in this method (Miller, in preparation), as both the high- 
frequency part of the PSD has a steep spectrum, P v oc v~ a where 
ol — 2.1 (McHardy et al. 2004) and also aliased high-frequency 
modes are suppressed by the sine 2 term, the effects of aliasing are 
minimal and we make no aliasing correction. 

The likelihood method may be extended to simultaneously 
measure the PSDs of two time series, and their cross-spectrum 
frequency-dependent coherence and time delays (Miller, in prepa- 
ration), and we show results from this analysis also. By incorpo- 
rating time delays as additional parameters in the likelihood model 
we can straightforwardly compute best-fit values and errors. 



4.2 Red noise leak and the PSD bandpower widths 

A well-known problem with any power-spectrum estimation 
method applied to finite time series is that modes on timescales 
longer than the window function are not correctly measured. In 
effect, long-timescale variations are interpreted instead as modes 
with periods on the scale of the window function, a phenomenon 
known as "red noise leak". If the variations don't have much power 
on long timescales this may not be important, but as NGC 405 1 has 
a low frequency PSD P(v) oc v~ a w here a ~ 1 at least down 
to v < 10~ 8 Hz jMcHardv et alll2004l) . red noise leak cannot be 
ignored. Miller (in preparation) ha s tested the effect of r ed noise 
leak using data simulated with the iMcHardv et alj ( 120040 PSD of 
NGC 405 1 and sampled with the window function provided by the 




log 10 v/Hz 



Figure 4. Comparison of the PSDs for NGC 405 1 derived from the three 
epochs: 2005 (solid curve), 2008 Nov 6-12 (dashed), 2008Nov 23-25 (dot- 
dashed), showing log 10 vP v against log 10 i^/Hz. Points with errors bars 
are plotted at the midpoint in log v of each bandpower, error bars show the 
Fisher matrix uncertainties, lines are drawn to connect the points in each 
epoch of observation. 



2008 Suzaku observations. As expected, it is found that the low- 
est frequency periods are measurably affected by red noise leak. 
In principle in our PSD estimator we can include low frequency 
modes by extrapolating the lowest bandpower down in frequency: 
however even this is not a good solution as in the data we are not 
averaging over a large number of low frequency modes, and fur- 
thermore the analysis becomes dependent on the model assumed 
for the extrapolation. 

A related issue is that modes close in frequency are statisti- 
cally correlated, because of the finite window function. We must 
therefore choose frequency widths for our bandpower estimates 
that are broader than the window function to avoid unpleasant cor- 
relations and degeneracies between neighboring bandpowers. At 
the lowest frequency bandpowers, the effect of red noise leak also 
introduces correlations and anti-correlations between neighbouring 
modes. The best solution is to choose bandpowers that are fairly 
broad. In the following analysis we choose bandpowers that are 
logarithmically spaced, with Alog 10 u = 0.4 or 0.5: at high fre- 
quencies this is adequate to ensure that they are independent. At 
low frequencies we bin together those bandpowers to make a single 
bandpower whose width is at least as large as the half-power full- 
width (HPFW) of the window function. We measure the HPFW by 
direct Fourier transform of the window function, and finding the 
the frequency lower than which half the summed power is found. 
The lowest-frequency bandpower in our analysis thus tends to be 
extremely broad in log v. 

Because of red noise leak, the PSD in this lowest bandpower 
should always be treated with caution in the subsequent analysis. 
One particular concern which arises is in the measurement of time 
lags between energy bands. A long period mode which has some 
small time shift between two energy bands will, in the low fre- 
quency limit, simply appear as having a shift in the mean value 
between the bands within the observed time window. The shift in 
the mean will be removed because we use the observed mean value 
to convert the measured time series into fractional fluctuations A, 
and hence the time lag will disappear. We return to this concern in 
sections l43H andl5T2l 



© 2009 RAS, MNRAS 000,m|T5l 



Suzaku X-ray spectrum of NGC 4051 7 




log 10 v/Hz log 10 v/Hz 



Figure 5. PSDs, plotted as vP{v) from the observations of NGC 4051 in 2005 (left) and 2008 (right). There are six energy bands: 0.3-0.94 (yellow), 0.94-1.2 
(magenta), 1.2-1.5 (cyan), 1.5-2.1 (green), 2.1-3.8 (blue), 3.8-10keV (red). The Fisher matrix uncertainties are shown as error bars. The PSD in the highest 
frequency bandpower in 2005 was consistent with zero in three of the energy bands, the plot shows the corresponding 2cr upper limits. 



4.3 Results 

4.3.1 The time series 

Time series were constructed from the data by binning and sam- 
pling at regular intervals. Any time bins when the source observa- 
tion covered less than half a bin were not used: the remainder were 
corrected by the "on-time" in each bin. As an illustration of the 
variability, Fig.[5]shows the full band, 0.3-10keV, time series sam- 
pled at 256 s (for the PSD determination in each band, sampling 
times of 64 s were used). The frequent and periodic gaps in the data 
are clearly visible. The variability is clear however, and even by 
eye it can be seen that it appears to have different characteristics 
in 2005 from the two observations in 2008. We quantify this in the 
next section. 



4.3.2 Non-stationarity of the variations 

The first Suzaku observation was obtained three years before the 
second two observations, when the source was in a substantially 
lower flux state. Even the second and third observations, taken 
two weeks apart, have different mean flux (sectionf2](. We there- 
fore start by comparing power spectra from the three observations 
to test whether they are consistent with stationary behaviour, or 
whether there are differences. Fig.|4] shows PSDs measured from 
each XIS dataset using the energy range 0.3-10keV, although be- 
cause of the steep source spectrum and the higher effective area 
in the soft band, the measured spectrum is dominated by the soft 
band. The time sampling used was 64 s, with a bandpower interval 
of A log 10 v — 0.4. Bandpower intervals smaller than the win- 
dow function full width (as defined in section l4~H were coadded as 
described above. The minimum frequencies differ for each dataset 
owing to the differing lengths of observation (TableQJ. 

The PSDs at low frequencies, v < 10~ J Hz, show significant 
variation owing to uncorrected red noise leak. Above this frequency 
however useful comparison may be made. The two observations 
from 2008 are closely consistent with each other, the third obser- 
vation having only a marginally significant excess of power at the 
highest frequency bandpower. The 2005 data appear significantly 
different however, with more power in the mid frequency range and 
a steeper spectrum to high frequencies. 



The mean levels used to normalise the fluctuations were dif- 
ferent in the three datasets of course. If we were to scale ob- 
servations 2 and 3 to the same mean, their amplitudes would 
differ. This implies either what has previo usly been suggested 
dUttlev & McHardvll200ll ; IUttlev et al.ll2005t) . that the fluctuations 
are intrinsically "fractional" in nature, or that some intrinsic fluc- 
tuations have been modulated by a secondary effect such as a 
varying multiplicative absorption. The 2005 observation, however, 
clearly has different shape. To quantify this, we measure the slope 
of the PSD at high frequencies by fitting a linear relation be- 
tween log P v and log v for v > 10 -4 Hz. Such a functional fit 
is a good fit for the 2005 and 2008 Nov 23-25 PSDs but fits less 
well the 2008 Nov 6- 12 PSD: as the aim of this exercise is to see 
whether there are statistically-significant differences between the 
PSDs of the three epochs, this simple parameterisation is sufficient. 
A minimum-x 2 -squared fitting method is adopted, where as well 
as taking account of the statistical uncertainty on each of the fit- 
ted PSD points, we also include the covariance between the points 
as determined from the Fisher information matrix. We find best-fit 
slopes a = 2.00±.08, 1.78±.03 and 1.67±.08, for each of 2005, 
2008 Nov 6- 1 2 and 2008 Nov 23-25 respectively, where P v oc u ~ a . 
Thus the two observations in 2008 have high-frequency slopes that 
are consistent with each other, but the 2005 PSD has a steeper high- 
frequency slope, that differs by 2.7cr from the 2008 Nov 6-12 slope 
and by 3a from the 2008 Nov 23-25 slope. 

We return in section|5]to discussion of these points in the full 
context of all the results, but for now we note that the PSD be- 
haviour of the source appears non-stationary, given the different 
shape PSDs between 2005 and 2008, and in the subsequent PSD 
analysis we treat the 2005 and 2008 data as separate time series. 



4.3.3 Energy-dependent PSDs 

To search for any dependence of PSD on photon energy, we divided 
the XIS photon events into six bandpasses covering the range 0.3- 
lOkeV. The likelihood method corrects for shot noise in the data, 
but as we shall be looking for differences between the PSDs of the 
six bandpasses, we have chosen energy boundaries such that each 
band has the same count rate. The source photon spectrum is steep, 
and the XIS effective area peaks in the soft band, so inevitably this 
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means the bands are not at all spaced evenly in energy. The energy 
bands were 0.3-0.94, 0.94-1.2, 1.2-1.5, 1.5-2.1, 2.1-3.8, 3. 8-10keV. 

The results are shown in Fi g.|5]for a time sampling of 64 s. As 
found by McH ardv et all {2004), and as expected from the longer 
timescale spectral variability in section[3] the harder bands have 
less overall variability. However, it is also clear that the hardest 
energy band, 3.8-10 keV, has a significantly flatter PSD than the 
softer bands, in fact having the highest power of any energy band at 
the highest frequency measured, despite having significantly lower 
power at lower frequencies, although given the error bars we cannot 
say that the hard band has significantly more power than the softer 
bands. T his result is in good ag reement with the model-fitting ap- 
proach of lMcHardv et alj J2004h applied to the 2001 XMM-Newton 
data, where evidence for flattening of the PSD with energy was also 
found. Overall, the shape and normalisation of the Suzaku PSD is 
in excellent agreement with that XMM-Newton PSD. There is some 
evidence that the lowe st frequencies dip belo w the best-fit relation 
P(u) oc v~ l found by McHardv et al. ( 20 0j), but the departure is 
consistent with the uncertainty in the McHardv et al .1 slope and in 
any case may be affected in our analysis by red noise leak. In the 
following analysis we ignore the PSD shape at v < 10 -5 Hz. 




T/s 

Figure 6. Time lags between the 0.3-0.94 keV band and the 6-10 keV (solid 
points and error bars) and 2. 1-3.8 keV (open points and dotted error bars) 
bands, as a function of time period T for the combined 2005-8 dataset. 
Vertical error bars show the 68 percent confidence regions. Horizontal bars 
show the range of time periods included in each bandpower. 



4.3.4 Time delays and coherence 

We can also search for evidence of time lags between energy bands, 
a phenomenon which has now been detected in a number of AGN, 
as discussed in the introduction. Again, we follow the method de- 
scribed by Miller (in preparation) where the cross power spectrum 
between two bands is calculated in terms of some amplitude and 
a time lag between the bands. The value of time lag derived cor- 
responds exactly to the time lag derived from the cross-correlation 
function in the case where every frequency bin is constrained to 
have the same lag. By working in the Fourier domain however 
we can also investigate frequency-dependent lags, effectively be- 
ing the lags obtained by filtering the cross-correlation function with 
a top-hat filter in frequency. The lags we derive by this method 
are directly compara ble to those derived by other authors (e.g. 
iMcHardv et al.ll2004h . but they have been arrived at by the maxi- 
mum likelihood method rather than periodogram or Monte Carlo 
methods. 

For this exercise we combine all datasets to maximise the 
signal-to-noise in the time lags. Even though the shape of the 
PSD varies between 2005 and 2008, the derived time lags show 
no significant evidence for epoch dependence: the 2008 data alone 
yield time lags very close to those from the combined dataset, 
with slightly larger uncertainties. The 2005 data alone yield time 
lags consis tent with the 2008 da ta but with large uncertainties. As 
found bv lMcHardv et al.l d2004l) . we find significant time lags be- 
tween energy bands, where the harder band lags the softer band 
and where the lag increases with energy difference. We define the 
energy bands similarly to section |4~3.3| Again, the bands have to 
be broad in energy in order to obtain an adequate count rate, but 
to best investigate the energy dependence of the time lags we re- 
lax the requirement to have equal count rates in all energy bands, 
and the hardest band from section |4~3.3| is split into two, 3.8-6 and 
6-10keV. The time sampling chosen was 256 s: this longer time 
sampling improves the count statistics in each time bin, and we 
shall see from the results that we are most interested in time lags 
on periods longer than 1 ks and with values for the lag measured 
in hundreds of seconds. Accurate characterisation of the measure- 
ment errors is important for the time delays, so here we supplement 
the standard Fisher matrix error measurements with direct measure- 



ment of the likelihood distribution, obtained by stepping through a 
range of time delays at each frequency, while maximising the like- 
lihood by allowing all other coefficients to vary. We define two con- 
fidence intervals by the regions 2A In £ ^ 1 and 2.7 respectively: 
as the likelihood distributions appear to be close to gaussian these 
correspond closely to the standard 68-percent and 90-percent con- 
fidence regions of a normal distribution. The 68-percent confidence 
regions are also close to the Fisher matrix errors, again as expected 
for gaussian uncertainties. Table|2] gives the Fisher matrix and 68- 
and 90-percent confidence regions for each frequency band. 

In these results, lags up to 970 s are detected with respect to 
the softest band, that increase with band energy and w ith the period 
of the fluctuations, as found by McHa rdv et al.l d2004l) . The longest 
period point is consistent with zero lag but with a very large un- 
certaint y. Even so, it does ap pear discrepant compared with the re- 
sults o f McHard v et al.l d2004l) and is worth further consideration. 
IMcHardv et alj found a lag of 3000 s on timescales ~ 70 ks be- 
tween the 0.1-0.5 keV and 5-10keV bands. The disagreement with 
the 90 percent confidence upper limit of 600 s (Table|2j could indi- 
cate a variation in source structure between 2001 and 2008: a more 
likely possibility, however, is that the long timescales have not been 
sufficiently sampled by the Suzaku observations to allow accurate 
measurement, in contr ast with the RXTE observations analysed by 
McH ardv et ail d2004l) that spanned 6 years. As discussed in sec- 
tion l4,2l red noise leak into the longest period bin may corrupt the 
measurement of the time delay. We therefore exclude this bin from 
further consideration. Time lags between the 0.3-0.94 keV band and 
the 2.1-3.8 and 6-10keV bands are shown in Fig.[6j plotted log- 
arithmically (other bands are not shown, for clarity, but may be 
inspected in Table[2} . 

We also investigated time lags between the PIN and XIS bands, 
however the PIN data were too noisy for any accurate time lag to be 
measured (although were consistent, within the large uncertainties, 
with the lags for the hardest band shown in Table[2j. 

Finally, we investigate the coherence, 1Z, between these 
two energy band s , wh ere coherence is defined following 
iNowak & VaugharJ (l996) as the ratio of the cross-spectrum am- 
plitude squared to the product of the individual PSD amplitudes for 
the two bands, as a function of frequency. It is defined to lie in the 
range ^ 1Z ^ 1. High coherence implies the different wavebands 
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Figure 7. Coherence TZ between the 0.3-0.94 keV and 6-10keV bands as a 
function of time period T for the combined 2005-2008 dataset (solid sym- 
bols). Vertical error bars show the Fisher matrix errors, horizontal bars show 
the range of time periods included. Also shown are the coherence between 
the same bands for the 2005 data alone (open symbols and dashed error 
bars) and for the 2008 data alone (cross symbols and dashed error bars). 
The highest time period in the 2005 dataset is not shown as the error bounds 
are too large to be useful. 



are highly correlated, low coherence implies they are not. Fig.[7] 
shows 1Z again as a function of frequency, showing only the Fisher 
matrix errors. The errors are large, so we have not directly mea- 
sured the likelihood distributions in this case (because we only plot 
the Fisher errors, the error bars are not constrained to lie within the 
range ^ 1Z ^ 1). Results from the highest time period in 2005 
are not shown as the uncertainty in 1Z is of order unity. In the com- 
bined dataset, 2005-8, and in the 2008 data alone, the coherence is 
high at all frequencies. However, investigating the 2005 low state 
data alone, the coherence drops for short period variations. 



4.4 Summary of key timing analysis results 

We also now summarise the key timing analysis results. 

(i) The PSDs demonstrate non-stationary behaviour between 
2005 and 2008. 

(ii) The PSD measured for the hard band is significantly lower 
in amplitude and with a flatter frequency-dependence than the soft 
band PSDs, as previously found by other authors. The difference 
between bands is more marked in 2005 than in 2008. 

(iii) We confirm the existence of both frequency-dependent lags 
of harder-energy bands with respect to softer-energy bands, up to 
970 ± 225 s, and the energy dependence of the lags. 

(iv) The coherence between soft and hard bands is high in the 
2008 data, but is low in the 2005 data for short period variations. 



5 DISCUSSION 

5.1 Comparison of timing results with previous observations 

The power spectrum and frequency-dependent time lags measured 
from the Suzaku observations of NGC 4051 are in remarkably 
close a greement with the measurements reported by McHardv et al. 
(2004), despite being based on data taken with different observato- 
ries (Suzaku v. XMM-Newton and RXTE) on widely separated dates 



Table 2. Time lags between the 0.3-0.94 keV band and bands at higher 
energy, as a function of time period T for the combined 2005-8 dataset. 
Columns give the range of periods in each band and the lag, followed by 
three estimates of the errors: Fisher matrix, direct evaluation 68% confi- 
dence region and direct evaluation 90% confidence region. The longest time 
period in each energy band is likely affected by red noise leak and is not 
used in the analysis. 



T range 
/ks /ks 



lag Fisher 68% 90% 
/s error /s error /s error /s 



6-10 keV 

2879 51.2 • 



28.8 
9.1 
2.9 
0.91 



16.2 
5.12 
1.62 



162000 
51.2 
16.2 
5.12 



0.512- 1.62 



-564 
970 
299 
66 
30 



±913 
±212 
±84 
±35 
±15 



+694 
-715 
±225 

±86 
±35 
±17 



+ 1163 
-1190 
+372 
-374 
±142 

±58 
±28 



3.8-6 keV 

2879 51.2- 162000 -872 ±711 

28.8 16.2-51.2 700 ±180 

9.1 5.12- 16.2 94 ±59 

2.9 1.62- 5.12 57 ±24 

0.91 0.512- 1.62 30 ±9 



+ 558 
-563 
±193 

±58 
±24 
±10 



+932 
-942 
±320 

±95 
±39 
±16 



2.1-3.8 keV 

2879 51.2 • 



28.8 
9.1 
2.9 
0.91 



16.2 
5.12 
1.62 



162000 
51.2 
16.2 
5.12 



0.512- 1.62 



-880 
429 
57 
28 
12 



±566 
±124 
±44 
±16 
±7 



+474 
-416 

±133 
±44 
±16 
±7 



+767 
-713 
±219 

±72 
±26 
±12 



1.5-2.1 keV 

2879 51.2 



28.8 
9.1 
2.9 
0.91 



162000 
16.2-51.2 
5.12- 16.2 
1.62-5.12 
0.512- 1.62 



-420 
373 
57 
7 
2 



±477 
±101 
±35 
±13 
±6 



+398 
-351 
±106 

±35 
±13 
±6 



+645 
-600 

±176 
±57 
±22 
±11 



1.2-1.5 keV 

2879 51.2 • 



28.8 
9.1 
2.9 
0.91 



16.2- 
5.12- 
1.62- 
0.512 



162000 
51.2 
16.2 
5.12 
- 1.62 



-444 
287 
48 
-1 
12 



±441 
±101 
±33 
±12 
±5 



+343 
-366 
±106 

±32 
±12 
±5 



+ 582 
-601 

±176 
±52 
±20 
±9 



0.94-1.2 keV 

2879 51.2 



28.8 
9.1 
2.9 
0.91 



162000 
16.2-51.2 
5.12- 16.2 
1.62-5.12 
0.512- 1.62 



-427 
216 
49 
10 
6 



±406 
±108 
±28 
±11 
±5 



+317 
-339 
±114 

±27 
±11 
±5 



+537 
-558 

±189 
±45 
±19 
±9 



(2005-8 v. 1996-2002) and using different analysis methods (max- 
imum likelihood versus Monte Carlo). The key results that are in 
common are: 

(i) NGC 405 1 is overall less variable at high energy than at low 
energy over most of the frequency range. 

(ii) The high-energy PSD is flatter so that by time periods about 
128 s the variability in each band is similar. 

(iii) The softer bands in particular show a sharp break in the 
power spectrum at ~ 10~ 4 Hz. 

(iv) Frequency-dependent time lags are found, with hard band 
lagging the soft band, with lag increasing with energy difference, 
and lag increasing with period. The maximum lag found in the 
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Suzaku analysis is 970 ± 22 5 s for periods of 10s of ks, similar 
to the maximum lag found by McHar dv et al.l d2004h on similar pe- 
riods of approximately 800 s. 

There are some differences: McHard v et al. N2oo4) find the las con- 
tinues to increase with period, which may be attributable to the 
much longer time baseline afforded by the 6 years of RXTE data 
they included; and they found a more significant decrease in co- 
herence with decreasing period which may be due partly to the 
difference in ana lysis methods a nd partly to the inclusion of low- 
state data in the M cHardv et alj analysis. By analysing the 2005 
low state separately from the 2008 higher state, we find significant 
differences in the PSDs and significantly lower coherence between 
bands in 2005. 



5.2 Reverberation time delays 

The frequency-dependent time delays in this and other AGN have 
previously been interpreted in terms of a model of intrinsic vari- 
ations hypothesised to explain those features, comprising pertur- 
bations propagating inwards on an accreti on disk where harder X - 
ray emission is emitted from smaller radii dArevalo & Uttlev| 2006). 
While such models are not ruled out by the analysis presented here 
(the results largely agree with previous analyses), they do not have 
any predicted relationship to the PCA spectral decomposition that 
we have seen in section[3] and although the earliest results indi- 
cated lags that varied linearly with time period, more recent results 
reveal a more complex dependence of lag on period that has been 
interpreted as indic ating the presence of multiple zones of fluctu- 
ating propagations (McHardy et al. 2007). However, there is a nat- 
ural explanation of the frequ ency-dependent time delays in terms 
of reverberation delays (e.g. |Peterson||l993l) which we can relate 
directly to the PCA spectral components. In section[3]we saw the 
existence of a hard component spectrum that dominates the low 
state, in common with many other AGN. If this component has 
a significant contribution from reflection, this may cause the ob- 
served time delays. In fact, reverberation models predict just the 
sort of frequency-dependent time delays that are observed here. We 
illustrate this with a simple model. 

Let us suppose that at any frequency, the continuum flux we 
observe comprises two components, one a component from a cen- 
tral source that is seen directly, plus a time-delayed component 
reflected from distant material. Consider one Fourier mode of the 
measured source variability with angular frequency uj, 



light per unit radius within a range ri r ^ r%, has an observed 
phase term given by 



a„,e 



b T e l " {t - r) dr, (4) 



where a u is the amplitude of that mode, r is the time delay of re- 
flected emission from some region, primed quantities denote the 
observed amplitude and time delay in the combined direct plus re- 
flected signal, b T is the transfer function defined over the range 
^ b T ^ T max and / is the total reflected fraction of light. We ex- 
pect ^ / < 1 unless the illuminating source is obscured. It can 
be seen that the phase of the combined signal, contained in the term 
e -»ujT ^ j s re [ ate( j to me F our i er transform of th e transfer function 
(see related results discussed by iPetersonll 19931) . A simple case is 
given by considering reflection from a thin spherical shell of radius 
r: the transfer function is then a uniform probability distribution 
over the range ^ r ^ r s where r s = 2r/c, and the imaginary 
part of equation[4] that is responsible for the observed time lags, is 
oscillatory. Smooth time lag functions may be produced by thick 
shells of reflection. A thick spherical shell, with uniform reflected 
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where we integrate over the thickness of the shell, defined by 
ri = 2ri / c and Irijc. Fig.[8]illustrates the predicted de- 

lays, t', as a function of the fluctuation time period T = 27r/oi 
for / = 0.8, Ti — 500 s and a range of values of r max . The 
frequency-dependence arises because high frequency modes are av- 
eraged out by the transfer function, and the remaining variability is 
thus more dominated by the illuminating source alone. Echoes of 
lower frequency modes can survive in the reflected signal, resulting 
in a measurable delay. The suppression of high-frequency modes in 
the reflected spectrum has an effect on the PSD, with the power in 
modes at uj > 1 /r max being suppressed by a factor (1 + /) ~ 2 with 
respect to those at lower frequencies. The coherence between the 
direct signal and the combined direct plus reflected signal tends to 
unity at low and high frequencies and varies as cos 2 (ci;r') at inter- 
mediate frequencies. 

One particular point to note is that the measured time lags are 
significantly less than the light travel time to the reflecting region, 
because the coaddition of primary and reflected light dilutes the 
time lag signature. In the limit to — > 0, the observed time lag tends 
to a maximum, 



/ ( r max Tl) 



(6) 



4(1 + /) (r max -ri)' 

as seen in Fig.[8]at high time periods for the models with low r max , 
so t' < r max always, and in estimating the distance of the rever- 
berating zone from the illuminating source we must allow for this 
dilution of the lag time. 

It would be interesting if a signature of a maximum time lag 
could be found in the data. In NGC4051 there is no evidence 
that su ch a maximum has be en found, either in this analysis or in 
that o f McH ardv etaD d2004h. However, a maximum lag ha s been 
found in Ark 564 dArevalo et alj200dlMcHardv et alj2 007). lend- 
ing support to the reverberation explanation. Another feature of the 
Ark 564 analysis is that at the highest frequencies there is some 
evidence for negative time lags (i.e. hard band leading soft band, 
McH ardv et~afl l2007) - in the reverberation model this may arise 
at high frequencies as a result of the time-lag oscillatory behaviour 
produced by a reflecting region that only partially covers the source. 
In effect, such apparent negative lags are an artefact created by the 
Fourier decomposition of the time lag signature. 

To measure r' in data, we need to be able to compare the com- 
bined, direct plus reflected, signal with the direct signal alone. In 
our measurements we find a lag between hard-band photons and 
soft-band photons. If we assume that in the soft-band we are seeing 
little reflected light, but that the hard band signal is a combination 
of direct and reflected light, we may compare Fig.[6]with Fig.U It 
can be seen that there is good qualitative agreement between the 
simple reverberation model and the observed time lags, although 
the time lag for the shortest periods is larger than expected. The 
profile of the frequency-dependence (the "lag spectrum") depends 
on the assumed transfer function: a fall-off in reflection with radius 
as 1 jr, in either spherical or flattened geometries, tends to produce 
a flatter relation between lag and period. As a simple example of 
a more complex transfer function, we consider a two-component 
reflection model, in which an inner reflecting zone has parameters 
/ = 0.3, ti — 200 s, r max = 800 s and a second zone extending to 
larger radii has parameters / = 0.8, n = 800 s, r max = 30000 s, 
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Figure 8. Predicted time lags, r', from the simple thick-shell reverbera- 
tion model as a function of fluctuation time period T for r max = 

10 3.5 S 

(dashed curve), f 4 s (solid curve) and 10 4 ' 5 s (dot-dash curve). The reflec- 
tion factor / = 0.8 for all curves. 




T/S 

Figure 9. Predicted time lags, t', from the two-component thick-shell re- 
verberation model, with parameter values as described in the text, as a func- 
tion of fluctuation time period. Also shown are the measured hard-band lag 
values from Fig.[6] Note that this model is illustrative only and has not been 
fitted to the measured values. 



where these parameter values have been chosen purely as an illus- 
tration and have not been obtained from fitting to the data (there 
are three parameters that primarily affect the shape of the lag spec- 
trum, and only four data points, so we do not consider fitting of this 
model to be particularly useful for this dataset). In this particular 
illustration, the total reflected light exceeds the light that is directly 
received by the observer, implying some loss of direct light by ei- 
ther absorption or scattering. The resulting predicted lag spectrum 
is shown in Fig.[9] together with the measured hard-band time lags 
from Fig.[6] Such additional complexity in the reverberation model 
is able to reproduce more complex lag spectra, leading to the hope 
that improved datasets may be able to better constrain the distri- 
bution of circumnuclear material. Although accurate measurement 
of the lag spectrum could allow the transfer function to be deter- 
mined, the present quality of data are insufficient to allow mean- 
ingful determination of the transfer function shape, and despite the 
lag excess on short timescales, overall the simple thick-shell model 
provides a reasonable description of the present dataset, and we 
now investigate the consequences of that simplest model. 

To place statistical limits on the parameters of the simple 
thick-shell reverberation model, we used the frequency-dependent 
lag model of equation|5] within the likelihood maximisation proce- 
dure to find best-fit values of / and r max , stepping through a grid 
of /, r max values and evaluating the best-fit likelihood for each 
grid parameter set. The value of the inner radius of reverberation 
is not well constrained by the data, and this was fixed at a value 
Ti = 500 s. Confidence levels were assigned by finding intervals 
in log(likelihood), A\ 2 = 2Aln£, of 2.3 and 4.61 correspond- 
ing to 68- and 90-percent confidence intervals. The best-fit values 
for the lags between the 6-10keV and 0.3-0.94 keV bands were 
/ = 0.8±g:Sg. log 10 r max = 4.0l°;J 5 (i.e. r max = 10 000l«°° s), 
where quoted errors are 68 percent confidence intervals, with a 90- 
percent confidence range in r max of 2 500 < r max < 20 000 s. The 
confidence intervals are shown in Fig.[l0] As r max = 2rifc, the 
reverberating region extends to a distance V2 = 5 000^ 2500 light- 
seconds in this model. We can see from consideration of Fig.|9]that 
the double shell model produces more reflection from smaller radii 
(in the example shown, up to r < 800 s) but with additional larger- 
scale reflection extending beyond the region where the current data 
can constrain its extent (r > 10 000 s). Irrespective of the detailed 
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Figure 10. 68- and 90-percent confidence intervals on the model param- 
eters log 10 T max and / for the simple single-shell reverberation model as 
described in the text, derived from fitting to the combined 2005-8 dataset. 
The cross marks the best-fit location. 

model, the continuous nature of the lag spectrum, lags monotoni- 
cally increasing with time period, implies there is thick-shell rever- 
beration over a wide range of radii. 

Overall, our conclusion is that reverberation naturally repro- 
duces the frequency-dependent lags. The reflection fraction / may 
be directly interpreted in terms of the global covering factor of the 
source, C' g : in the case where our line of sight to the source is un- 
obscured, / = 0.8 implies the remaining sight lines are 80 percent 
covered by material if it is 100 percent reflective, when integrated 
through the reflecting zone, with higher covering factors for lower 
integrated reflection albedo. However, the inferred global covering 
factor is smaller if our line of sight is partially obscured, as this 
means we could have overestimated /. In particular, if the cover- 
ing fraction of our line of sight is the same as the global cover- 
ing factor, and if that material transmits no light but has an albedo 
A, we expect / ~ AC g /(l — C g ) from which we would deduce 
C g = f/(A + /). Thus for / = 0.8 we deduce C g > 0.44 for 
A < 1. 

The reverberation time lags we deduce in either the single- 
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Figure 11. The variation in reverberation model reflection factor, /, (the 
ratio of reflected to directly-received flux) with the median energy E of 
each band. 

shell or more complex models place a significant component of 
reflecting material well within the broad-line region. In the case of 
the single shell model, it extends to distances r-z — 1.5 x 10 14 cm, 
or 600 x g (where r 9 = GMbh/c 2 ) with a 90 percent confidence 
range of 150-1200 r 9 for the nominal black hole mass of Mbh = 
1.7 x 10 6 Mq. More complex models may have a component of re- 
flection extending to larger radii, but this must be counter-balanced 
by enhanced reflection also at smaller radii. The most likely reflec- 
tion configuration is discussed further in section lBTol 

5.3 The reflection spectrum 

McH ardv et all d2004l) found that the time lass increase with the 
difference in energy between the bands being cross-correlated: 
harder bands are more lagged with respect to softer bands. This is 
confirmed by the analysis here. In the reverberation model, this ob- 
servation implies that the fraction / of reflected light increases with 
increasing energy. We can quantify this by fitting the simple "single 
thick shell" reverberation model to the time lags measured in vari- 
ous bands, and for a fixed model value of r max we can deduce / as 
a function of energy band, again assuming that the softest band has 
a negligible reflection component. This is shown in Fig.QT]for the 
same energy bands used above, for the best-fit value r max = 10 4 s. 
The median photon energy in each band is plotted. 

It is striking how similar the energy dependence looks to the 
spectrum of the offset component in the PCA analysis, which we 
have already interpreted as arising from a combination of absorp- 
tion and reflection. We thus have a natural explanation, not only 
of the frequency dependence of the lags, but of their energy de- 
pendence also. The fall-off to low energy is a direct consequence 
of scattering in a medium whose opacity increases with decreasing 
energy. It is a generic feature of mode l s of optically-thick reflec- 
tion (e.g. lMagdziarz & ZdziarskHll995URoss & Fabianll2005l) . but 
is expected to be seen also in the spectrum of scattered light in a 
less optically thick medium. A full calculation of Compton scatter- 
ing in an extended, ionised, possibly clumpy, medium is required 
to make a full prediction of the expected spectrum. Of particular 
interest is the strength, line energy and profile of Fe Kq emission 
that would be expected. We saw in section POl that there exists a 
broadened component of Fe Ka on eigenvector one with equiva- 
lent width 66I37 eV. If the reflected light accounts for a fraction 
//(l + /) of the total observed light, for / = O.St'.H we find 



the equivalent width against the reflected component alon^f] would 
be about 150±gs 3 eV. This is significantly lower than would be ex- 
pected for reflection from neutral material, where equiva lent widths 
of order 1 keV are expected (e.g. lGeorge & Fabianll99lt) . but in the 
environment close to the primary X-ray source we expect gas to be 
highly ionised, and such equivalent widths are produced in models 
of reflection in ionised mat erial, for values of ion isation parame- 
ter ^ > 100 erg cms -1 (e.g. iRoss & Fabianll2005h . Two particular 
effects that would need to be taken into account in modelling the 
expected line emission are the effects of Compton scattering in a 
hot medium, which could significantly broaden the line, and the ef- 
fects of resonant absorption and resonant line destruction in m ate- 
rial where Fe is io nised to higher than Fe XVII jMatt et alj 19971 . see 
also the review bv lTurner & Millerll2009l and references therein). 

In practice, the reflection spectrum and the dependence of time 
lag on frequency would be expected to have a more complex inter- 
relationship than discussed in this simple model. Material close 
to the ionising source would be expected to have a higher ioni- 
sation, and hence be more reflecting at lower energies, than mate- 
rial further away. Long period fluctuations would thus have reflec- 
tion spectra characterised by lower ionisation than the reflection 
associated with short period fluctuations. The structure and geom- 
etry of the reflecting region could also have a significant effect on 
the reflected spectrum: for example, a clumpy two-phase medium 
could simultaneously produce a "grey" reflection spectrum from 
fully ionised material and a hard reflection spectrum from lower- 
ionisation clumps. Or a highly ionised reflecting region might itself 
be surrounded by lower ionisation absorbing material. Such models 
could produce a wide range of spectrum hardness and emission line 
equivalent widths. We postpone investigation of more sophisticated 
models to future work. 

5.4 Coherence and the amplitude of the PSD 

The most noticeable feature of the comparison of the 2008 PSD 
between hard and soft bands is that the high energy bands show 
less fractional variability than the soft bands (i.e. the amplitude of 
the PSD is significantly higher in the soft bands compared with 
the hard band). The coherence in the overall dataset is high, how- 
ever, indicating that the same fluctuations are present in soft and 
hard bands but that their amplitude in the hard band is suppressed. 
Strangely, though, at the highest frequencies the soft and hard bands 
have comparable PSDs, with the hard band variations possibly be- 
ing even higher than the soft band, which does not fit easily with 
the simplest interpretations. 

First, consider the expectation of the reverberation model. In 
the combined direct plus reflected signal the power spectrum is 
expected to be suppressed by a factor (1 + f) 2 at frequencies 
uj > l/r max , compared with the direct signal alone (section [5~2l . 

2 One caveat is that we need to be careful in interpreting the PCA com- 
ponents in this context: they were derived from spectral variability on 
timescales longer than 20 ks, but in the reverberation signal we are studying 
variations on much shorter timescales, and we do not know that the PCA 
spectral decomposition applies on these timescales. If the long-period vari- 
ations that shape the PCA are caused by partial-covering variation in the 
inner regions, some component of the reflection spectrum from the rever- 
berating region could appear in the offset spectrum rather than being en- 
tirely contained within eigenvector one, and thus the true equivalent width 
of Fe Kq on the reflection component could be higher than estimated here. 
Unfortunately we lack sufficient signal-to-noise to apply the PCA method 
to faster variations while retaining the spectral information. 
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Thus if the soft band PSD is dominated by intrinsic variations, and 
the hard band PSD is a reverberation-modified version of that, then 
we expect the 2008 hard band PSD to be lower than the 2008 
soft band PSD by a factor about (1.8±;J|) a ~ 3.2±;|. This is 
consistent with the difference in PSD amplitudes at frequencies 
10" 5 < v < 10" 3 Hz (Fig.g} for r max > 10 000 s. 

However, this expected PSD amplitude difference is not ob- 
served at the highest frequencies, v > 10~ 3 Hz, in 2008. Here, the 
hard band had at least as much variance as the softer band. It has 
previously been suggested that the flatter hard-band PSD is a re- 
sult of the mechani sm that produces the intrinsic fluctuations (e.g. 
lArevalo et ai1l2006h . which we have not addressed here. If we be- 
lieve that the fluctuations we observe are intrinsic, the flatter PSD 
could be a consequence of that process. As an example of a possible 
effect that could cause it, the highest frequency modes are likely not 
independent gaussian fluctuations, but are more likely the high fre- 
quency modes associated with sharp spikes of emission, so in this 
picture soft-band spikes might last longer than hard-band spikes. 
Without a more specific model for the emission process it is diffi- 
cult to make further progress on understanding the PSD. A future 
statistical analysis of non-gaussianity in the time series may shed 
more light on the high-frequency behaviour. 

The coherence of the hard and soft band signals may also 
aid our understanding. In 2008 the coherence was high at all fre- 
quencies, with a minimum coherence TZ ~ 0.5 at intermediate fre- 
quencies, rising to higher coherence at high, and probably at low, 
lo (Fig.[7j. Reverberation is qualitatively expected to show simi- 
lar behaviour (section [5T2l (. but as the maximum observed value 
of ujt' ~ 0.21 (TableO, the minimum coherence expected from 
reverberation alone is TZ = cos 2 (ut') > 0.95. Reverberation 
might make some contribution to the dip in TZ at intermediate fre- 
quencies, but overall the observed coherence values in the range 
0.5 < TZ < 0.9 are more likely due either to imperfect intrinsic cor- 
relation between source fluctuations in the soft and hard bands or to 
some weak additional component of variability such as might arise 
from variable obscuration. Such effects cannot dominate, however, 
otherwise they would destroy the reverberation time lag signatures. 

Turning to the 2005 low-state data, the PSDs have a different 
shape from the PSDs in 2008 (section |4.3.3b and show both en- 
hanced soft band variability on longer timescales and reduced co- 
herence between the soft and hard bands. These two observations 
are indicative of an additional source of variability affecting the soft 
band but less-so the hard band in 2005. Given the spectroscopic ev- 
idence that the central source was obscured in 2005 (section[3]l, a 
strong candidate for that additional variability is varying absorp- 
tion, either changes in column density and ionisation, or variations 
in covering fraction. A varying covering fraction of a clumpy ab- 
sorber would match the longer timescale spectral variability uncov- 
ered by the PCA. In the hard band, some variations arising from the 
absorber changes would be expected, but greatly reduced in am- 
plitude because of the lower opacity at higher energies. As these 
variations would be uncorrelated with any intrinsic variations, the 
soft and hard bands would become more decorrelated, as observed 
(Fig.0. 

5.5 The absence of inner-disk reflection signatures 

In Section[3]we saw how there exists a hard component of emission 
which appears to be quasi-steady, and essentially represents the low 
spectral state of the source. It is this component which, in mean 
spectra, leads to the existence of a "red wing" at energies below the 
6.4keVFeKa line, and which has previously been widely inter- 



preted as being highly redshifted reflection from the accretion disk 
close to or even inside the last stable orbit of the black hole. The 
general lack of variability of this component is problematic how- 
ever, as its amplitude should vary with the ampli tude of the illumi- 
nating continuum. To counter this expectation, Fabian & Vaughar] 



( 120031) . iMiniutti et all d2003h and IMiniutti & Fabianl fc004h have 
proposed a model in which the continuum source is compact, and 
lies close to the black hole. The source is hypothesised to move 
towards and away from the accretion disk, as it does so light rays 
follow geodesies in the space time, and with careful placement it 
can be arranged that light reflected from the disk and reaching the 
distant observer can have an approximately invariant flux, while the 
flux from the source itself is low when the source is near the disk 
and high when it is further away. 

The large time lags found here must however come from more 
distant material, and even though in NGC4051 there is a compo- 
nent of hard excess that varies with the primary illumination, the re- 
flecting material must be distant, not from the inner disk. We cannot 
rule out the existence of additional inner disk reflection, with essen- 
tially zero time lag, but given that we already need distant (lagged) 
reflection there is no justification for invoking an additional inner 
disk reflection component. Furthermore, the hard spectral compo- 
nent can be modelled quite satisfactorily as being caused by re- 
flection and/or absorption from in tervening material (e.g. in MCG- 
6-30-15, [Miller et al J [20081 120091 see also iTerashima et alj|200^ 
and Lobban et al. (in preparation) for discussion of NGC4051). 
In the reverberation signal we find evidence for a high global cov- 
ering factor of reflecting components, so all spectral features are 
accounted for. 



5.6 Towards a complete picture of spectral variability 

In this final discussion section, we attempt to collate the various 
observational results and interpret them in terms of the reverbera- 
tion picture that has emerged in the previous sections. The aim is 
to build a complete model of the source spectral variability that is 
consistent between the spectral and timing analyses. 

First, the difference between the 2005 and 2008 spectral states 
may be explained as being caused by obscuring material block- 
ing a direct view of the central source, so that the flux we see is 
dominated by reflection in 2005. This explains the strong invari- 
ant 6.4keVFeKa line as well as the spectrum continuum shape 
if the line-emitting material continued to receive an invariant ion- 
ising flux: that is, if most lines of sight other than ours remained 
open. Low, reflection-domi nated states of NGC4 51 have pre- 
viousl y be en identified by Guainazzi et alj dl998l) , IPounds et al.l 
(2004) and Terashima et al. (2009). The obscuration idea was pro- 
posed bv lPounds et alj (2004) before the Suzaku da ta became avail- 
able, a nd i s confirmed by the late r observations. iGuainazzi et al.l 
dl998h and ITerashima et al J fe009l) also proposed an invariant re- 
flection component, although their suggestion that its constant na- 
ture is due to light travel time effects is not required if most of the 
transition between low and high flux states is caused by obscur- 
ing material along our line of sight to the central source. Because 
we did not have a direct view of the central source in 2005, any 
rapid variations would have been erased by reflection travel time 
delays, causing the steeper cut-off to high frequencies in the 2005 
PSD compared with the 2008 PSD (Figsg]&[5j. We do not have 
any good constraint on the location of the constant reflection com- 
ponent, however. The low coherence between soft and hard bands 
on timescales < 5000 s and the large offset between the PSD am- 
plitudes for the soft and hard bands in 2005 implies the presence of 
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a significant additional source of variability in the soft band, which 
we suggest is variable obscuration, most likely changes in covering 
factor, during the 2005 observation. 

Moving to the 2008 observations, the hard excess increased 
with the general brightening of the source. There are two possi- 
bilities. If the variable component of hard excess arises from re- 
flection, but if the 2005 state was indeed due to obscuration, the 
postulated variable reflection component must also have been sup- 
pressed in 2005, i.e. it was also partially obscured along our sight- 
line. This would be co nsistent with the partial-covering fitting by 
iTerashima etafl d2009h . There is evidence for a moderately broad 
component of Fe Ka that varied with the variable continuum (sec- 
tions l3.2l and l5.3l >. The variable hard excess and associated FeKa 
component would need to come from a region closer in than the re- 
gion responsible for the steady, narrower Fe Ka, so that it could be 
obscured in this way. Alternatively, in 2008 we could still have been 
looking through some partially-covering absorbing ma terial that al- 
lowed transmission at high energy (as in PDS 456, iReeves et al.l 
2009) and which was suppressed in 2005 because of a higher col- 
umn density. In either case, we require the hard excess to be par- 
tially obscured in 2005 and possibly in 2008 also. Both possibilities 
require obscuring and reflecting material to exist over a wide range 
of radii. Since 1998, the Fe Kg streng th has been constant to within 
a factor two dTerashima et alj|2009h : as the emitting region may 
well be substantially smaller than 10 light years, the implication is 
that the intrinsic illuminating continuum likewise has been constant 
to within this factor, and that the larger amplitude continuum vari- 
ability on long timescales is due to variations in obscuration along 
our line of sight. 

The time lags and their frequency- and energy-dependence in- 
dicate reverberation from reflecting material. The fraction of re- 
flected light must be high, / ~ 0.8 at E > 6 keV, indicating a high 
global covering factor for this gas. Model fits of a thick reflecting 
shell imply the reverberating zone is within a few light hours of 
the central source, extending over distances up to 5 000^ 2500 light- 
seconds from the illuminating source. The spectrum of reflected 
light required to match the energy dependence of the lags is con- 
sistent with having the hard spectrum generically expected from 
reflected emission from ionised material. Obscuring and reflecting 
material may, and probably does, exist at larger distances, but we do 
not expect to see signatures of longer time lags in the timing anal- 
ysis presented here, as longer timescales would need to be probed. 
We note that reverberation signals are not created by variations in 
obscuration along the line of sight alone, because most sight lines 
would be unaffected. Thus the primary fluctuations whose echoes 
we see in the reflected light must either be intrinsic or, if due to 
varying obscuration, must be such that the majority of sight lines 
are affected. 



The implications so far then are that long-timescale varia- 
tions may be caused by changes in obscuration along the line of 
sight (likely to be partial covering changes) but that variations on 
timescales of days or less are most likely caused by an intrinsic 
mechanism. Because the coherence between energy bands is high 
over the observed timescales of days and less in 2008, and because 
we expect varying obscuration to be uncorrelated with intrinsic 
variations, there is not much room for both varying obscuration and 
intrinsic variability on these timescales in the more common high 
states of the source, such as in 2008. 



6 CONCLUSIONS 

We have considered both the spectral variability and the broad-band 
variability power spectrum of three Suzaku X-ray observations of 
NGC4051 spanning three years. We find and confirm many fea- 
tures from the spectral analysis and the power spectrum analysis 
that have been found by previous authors using independent data 
and methods, both in this AGN and others. We have related the 
rapid variability features to the longer-timescale spectral proper- 
ties, and we have shown that a self-consistent picture between these 
regimes is possible. 

The interpretation of the long-t erm spectral v ar iabilit y 
of NGC4051, pre viou sly noted by Guainazzi et al] dl998h . 
Pounds ct iD d2004b and ITerashima et alj d2009T) . as being due to 
the primary X-ray source disappearing from view, leaving a con- 
stant, hard, reflection-dominated component, is confirmed here by 
the principal components analysis, and we reiterate the suggestion 
that this is due to variable partial covering obscuring the central 
source during observed low states, explaining the almost constant 
Fe Ka narrow line flux. The PSD differs in the 2005 low state from 
that in the 2008 higher-state, showing features consistent with this 
explanation. 

We have detected the previously-known time lags between 
hard and soft bands in the new Suzaku X-ray light curve of 
NGC4051. Although previously hypothesised to be in some way 
related to the mechanism that produces X-ray emission, here we 
conclude that these time delays, and by inference those in other 
AGN also, are well described by reverberation from material up 
to a few light hours from the central source. The important signa- 
tures are a strong dependence on the frequency of variations being 
measured and on photon energy, both of which have a natural ex- 
planation in the reverberation model. The frequency-dependence 
is a consequence of observing superimposed direct and reflected 
components, and the measured lag-frequency relation is consistent 
with reverberation from a thick shell of material, although the pre- 
cise geometry is not yet well constrained by the data. The measured 
energy-dependence implies a hard reflection spectrum that is con- 
sistent with that expected from ionised material and which is also 
similar to the hard component that is seen in principal components 
analysis of the spectral variations. To obtain sufficient reflection, 
the covering factor of reflecting material must be high, greater than 
about 0.44, as recently inferr ed from spectral modelling of other 
AGN (e.g. lMilleretal]|2009h . The existence of a maximum ob- 
served time lag at low frequencies and possib le negative time lags 
at high frequencies, as observed in Ark 564 dArevalo et alj|200q ; 
iMcHardv et alj|2007l) are expected in the reverberation model. 

The identification of the signal as being reverberation con- 
firms the idea that NGC 405 1 and similar AGN are surrounded by 
substantial covering factors of optically-thick, or nearly so, mate- 
rial within a few light-hours of the central source (far inside both 
the optical broad-line region and the dust sublimation radius). Al- 
though this AGN has a spectral hard-band excess very similar to 
that seen in other type I AGN, there is no requirement for any 
relativistically-redshifted reflection from the inner accretion disk, 
which would not generate the observed time lags. The reverberat- 
ing material extends over distances up to a few hundred gravita- 
tional radii from the illuminating source. Given the spectral and 
timing similarities with other type I AGN, we suggest that this 
source structure is a common feature of the type I population. 
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